Solving the heat-flow problem with transient relativistic fluid dynamics 
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Israel-Stewart theory is a causal, stable formulation of relativistic dissipative fluid dynamics. This 
theory has been shown to give a decent description of the dynamical behavior of a relativistic fluid in 
cases where shear stress becomes important. In principle, it should also be applicable to situations 
where heat flow becomes important. However, it has been shown that there are cases where Israel- 
■ Stewart theory cannot reproduce phenomena associated with heat flow. In this paper, we derive a 

relativistic dissipative fluid-dynamical theory from kinetic theory which provides a good description 
of all dissipative phenomena, including heat flow. We explicitly demonstrate this by comparing this 
. theory with numerical solutions of the relativistic Boltzmann equation. 

PACS numbers: 12.38.Mh, 25.75.-q, 11.25.Tq 

O. 

(T) ; I. INTRODUCTION 

The derivation of a consistent theory of relativistic fluid dynamics has been a challenge for some time. The difficulty 
resides in the parabolic nature of Navier-Stokes theory, which allows signals to propagate with infinite speed. While in 
non-relativistic theories this unphysical feature can be dismissed, in relativistic systems it leads to unstable equations 
of motion and is simply unacceptable |l| . A necessary (but not sufficient) condition is that the equations of motion 
are hyperbolic. 

For relativistic dilute gases, the Boltzmann equation provides a reasonable description of the underlying microscopic 
theory and can serve as a starting point to investigate the microscopic foundations of fluid dynamics. For the sake 
of simplicity, in this paper we consider a single-component system of massless and classical particles interacting only 
via elastic two-body collisions with a constant cross section a. The relativistic Boltzmann equation for this type of 
system reads 
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k^dj k = i J dK'dPdP' Wuu^pp, (/p/p, - / k /k0 , (1) 
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where /k is the single-particle distribution function, Wkk/-s-pp/ is the Lorentz-invariant transition rate, and dK = 

ge? 3 k/ [(27r) 3 fc°] is the Lorentz-invariant momentum-space volume, with g being the number of internal degrees of 
. . . freedom. The right-hand side of the Boltzmann equation describes the change of the single-particle distribution 

function due to elastic collisions of two particles with incoming momenta k and fc', and outgoing momenta p and p' . 
k> , The conservation laws follow from the first two moments of the Boltzmann equation as a consequence of the 

conservation of particle number and energy-momentum in individual collisions, d^N^ = 0, d^T^ — 0. The particle 
d \ four-current, iV M , and the energy-momentum tensor, T^ v , can be cast in the usual form, 

N^ = nu^ + n^, T» v = eu^u" - P Q /SJ 1V + -k^ , (2) 

where n is the particle number density, is the particle diffusion four-current, e is the energy density, it M is the 
fluid four- velocity defined in the Landau frame @, i.e., u v T^ w = eu^, Pq is the thermodynamic pressure, and 7r' Il/ is 
the shear-stress tensor. Since we consider a massless gas, the bulk viscous pressure is always zero. The additional 
equations needed in order to solve the conservation laws, i.e., the equations of motion for and tt^ , must be derived 
by properly matching fluid dynamics to the relativistic Boltzmann equation. 

In the framework of relativistic kinetic theory, Israel and Stewart were among the first to derive a relativistic theory 
of fluid dynamics that is causal and stable [3|. In the Israel-Stewart (IS) formulation, the single-particle distribution 
function is expanded around its local equilibrium value in terms of a series of Lorentz tensors formed of particle 

four-momentum i.e., 1, k^k v ', In order to derive fluid dynamics, this series is truncated at second order 

in momentum, i.e., one only keeps the tensors 1, fc M , and fc^fc", the so-called 14-moment approximation. The 14 
coefficients of this truncated expansion are uniquely matched to the 14 components of and T^ v . Finally, the 
equations of motion for and ~n lLV are obtained by substituting the truncated moment expansion into the second 
moment of the Boltzmann equation. 



2 



In the past five years, IS theory has been widely applied to ultrarelativistic heavy-ion collisions in order to describe 
the time evolution of the quark-gluon plasma (QGP) and the freeze-out of the hadron resonance gas appearing in 
the late stages of the collision. However, in heavy-ion collisions extreme conditions occur which question the validity 
of fluid dynamics. The QGP created at the Relativistic Heavy Ion Collider (RHIC) and, recently, at the Large 
Hadron Collider (LHC) is not only the fluid with the smallest space-time extension (~ 10 fm) ever created in nature 
but also the one where the space-time gradients of the fluid fields, for instance energy density e, are the largest 
(<~ Ic^el/e ~ 1/fm) ever encountered. On the other hand, Israel and Stewart's derivation lacks a small parameter, 
such as the Knudsen number, with which one can do power counting and systematically improve the approximation 
to describe higher-order gradients. Thus, the applicability of IS theory to the extreme conditions reached in heavy-ion 
collisions is, at the very least, not clear. 

One way to investigate the applicability of IS theory is to compare the solutions of this theory with numerical 
solutions of the Boltzmann equation So far, it was confirmed that, at least for some special problems, the IS 

equations Q are not in good agreement with the numerical solution of the Boltzmann equation. In these cases, it 
was shown by some of the present authors [f| that the IS formalism seems to be unable to describe heat flow, even 
when the Knudsen number is very small. In this paper, we would like to explain the reason for this discrepancy and 
propose a solution which works at least if the Knudsen is not too large. 

Recently, a systematic derivation of fluid dynamics from the Boltzmann equation was introduced in Ref . 9;] . The 
main difference between IS theory and the theory derived in Ref. Q is that the latter does not truncate the moment 
expansion of the single-particle distribution function. Instead, dynamical equations for all its moments are considered 
and solved by separating the slowest microscopic time scale from the faster ones. Then, the resulting fluid-dynamical 
equations are truncated according to a systematic power-counting scheme using the inverse Reynolds number R _1 ~ 
\n^\ jn ~ / Pq and the Knudsen number Kn = A m f p /L, with A m f p being the mean free-path and L a characteristic 
macroscopic distance scale, e.g. L -1 ~ d^u^ . The values of the transport coefficients of fluid dynamics are obtained 
by resumming the contributions from all moments of the single-particle distribution function, similar to what happens 
in the Chapman-Enskog expansion [Toj . 

In this paper, we review the basic aspects of this resummed transient relativistic fluid-dynamical theory (RTRFD). 
We explicitly derive the equations of motion of RTRFD including terms up to second order in the Knudsen number. 
We then demonstrate that this method is also able to handle problems with strong initial gradients in pressure or 
particle number density. This resolves the previously observed differences between the solution of IS theory and of 
the Boltzmann equation Q mentioned above. We conclude that these differences were caused by the uncontrolled 
truncation procedure of the expansion of the single-particle distribution function in terms of Lorentz tensors in four- 
momentum as employed in IS theory. 

This paper is organized as follows. In Sec. [Hi we review the basic aspects of RTRFD. One important point is that 
the fluid-dynamical equations of motion, as derived in Ref. [t| , become parabolic once second-order terms in Knudsen 
number are included. In Sec. IIII1 we explain how to obtain hyperbolic equations of motion in the RTRFD formalism 
up to second order in the Knudsen number. In Sec. IIV1 we compare solutions obtained within RTRFD for various 
levels of approximation and within IS theory with numerical solutions of the Boltzmann equation computed using 
BAMPS [111 ] . In Sec. |V]we discuss our results and draw conclusions. The Appendix contains intermediate steps of 
our calculations. We use natural units h = c = ks = !■ The metric tensor is g^ v = diag (+, — , — , — ). 



II. REVIEW OF RESUMMED TRANSIENT RELATIVISTIC FLUID DYNAMICS 

In RTRFD |9J, /k is expanded in terms of an orthonormal and complete basis in momentum space. The expansion 
basis contains two basic ingredients: The first are the irreducible tensors, 1, k^K k^ 1 k^ 2 \ . . . , k^ 1 ■ ■ ■ k^ m \ which 
form a complete and orthogonal set, analogously to the spherical harmonics 0, 

mm. Here > 

we use the notation 

^(w— = Afi*'"{ff A" 1 '""* , with A^.'.'.'y™. The latter quantities are projectors onto the subspaces orthogonal to u M . 
Their definition is explicitly given in Refs. Him- Except for m = 1, where A£ = g^ — u^u v , they are traceless. E.g., 
for m = 2, A^ = (A£A^ + A^A^)/2 - A^A^/3. Note that the expansion of / k in IS theory is not in terms of the 

irreducible tensors k^ 1 ■ ■ ■ k^ m ' , but in terms of the tensors k^ 1 ■ ■ ■ k^ m which are complete but neither irreducible 
nor orthogonal. 

The second ingredient are orthogonal polynomials in E^ = u^k^, pj^ — X}"=o ^or details in constructing 

the polynomials, see Ref. 0. Then, /k is expanded as 

A = /ok + / ok £E ntl pjr ,w ■ ■ ■ k^ , (3) 

£=0 n=0 

where /ok = exp (a?o — floE^) is the local equilibrium distribution function, with ao = ji/T being the ratio of chemical 
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potential to temperature and (3 = 1/T the inverse temperature. We further introduced the energy-dependent 

coefficients = (W^/£\) J2m=n "mn-ffki with a normalization constant and the irreducible moments of 
Sfk = fk — /ok, 

■^-M ee [ dK ££ k {fil ■ ■ ■ fc M 6f k . (4) 



J 

Some of the irreducible moments are related to the fields in Eq. @: — p^ and 7r Ml/ = p^ v . The values of ao an d 
(3q are defined by the matching conditions, n ee (£k)cb e = (^k)o' wnere (' ' 4 )o = / dK (• • • ) /ok- The matching 
conditions and the definition of according with the Landau picture imply that the following moments should vanish: 
Pi = 92 = Pi = 0. 

The equations of motion for p% and p^ v together with their respective transport coefficients were derived in Ref. [9| . 
Since we are investigating the massless limit, the scalar moments p r will not play a dominant role (they contribute 
mainly to the bulk viscous pressure) and we set them to zero. We also neglect all irreducible moments with tensor 
rank higher than 2, since they are traditionally not considered in fluid dynamics. The role of such moments will be 
investigated in a future work. Then, the equations for the first- and second-rank tensors pj? and p% v read 



p^ + £ A (2 „Vr = 2a( 2 >^ - 1 (2r + 5) p^ + 2p^> + f V^jx - \ (r + 5) p%^ T -±± p^6, 



n=0 



(5) 



where "' w> = A^y/.^u^d^-^, I» = V^ao, = u v \ and V = 9]. We also defined the thermody- 
namic integrals 

I ng (a ,A>) = {2q l m J dKEr 2q {-A^kckp) Vok . (6) 

The coefficients 

= ^-L_ j dKdK'dPdP'Wkk^ppJokfo^fopfop^lT 1 ^ ■■■k v ^ 

x (H&kfr ■■■k vt) + U<$ n k{ Vl ■■■k' vi) - H p e l P{vi ■■■p vt) - n% P ' {vl ■■■p' Vl) ) (7) 

contain all the information of the microscopic theory, while a-r are complicated functions of /3o and ao 0| . 

We now identify the microscopic time scales that dominate the long-time dynamics of the Boltzmann equation. 
This can be achieved by finding the normal modes of Eqs. ([5]), i.e., by diagonalizing the matrices A^ . Then, one 
obtains the following (in the linear regime, decoupled) set of equations of motion, 

xy> + xi 1] x? = pWi» + Xf, 
x^ v) + x^xr = Pl 2) ^ + K v > (8) 

where xf 1 are the eigenvalues of A^ and Xjf 1 '"^ 1 = Ylf=o (■^ _1 )rri Prn W are ^ e e ig enmo des of the linearized 

Boltzmann equation, with A^fl^ — diag(%o^> • • • > x¥\ ■ - •)• The terms Xjf and X^ v represent nonlinear 

terms and terms containing derivatives of the eigenmodes. The coefficients /3r are complicated functions of ao and 
Pa 01- Without loss of generality, we order the eigenmodes of the linearized Boltzmann equation, X^ 1 according 

to increasing Xr'\ i-e., in such a way that Xr^ < Xr+n ^ 

It is clear that when the terms X^ and Xfin Eqs. ((HJ are small, each eigenmode relaxes independently to its 
asymptotic solution on time scales given by 1/xr '■ We shall refer to the solution at asymptotically long times as 
Navier-Stokes value. The slowest varying eigenmodes of the Boltzmann equation are those for r = 0. At very long 
times, only the transient dynamics of these eigenmodes has to be resolved, while all quickly varying modes are assumed 
to have already relaxed to their corresponding Navier-Stokes values. In practice, this means we have to solve Eqs. © 
for r = dynamically, 

x^+x^x^pWiv + xZ, 

XW+^XF^PP^ + XP, (9) 
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while all other eigenmodes [r > 1) are approximated to be very close to their asymptotic (i.e., Navier-Stokes) values, 



Xtf 



xl 1 



■I" 



X£ v 



Xr 



r > 1 



(10) 



where the dots indicate deviations from the asymptotic solution which are at least of order 0(Kn 2 ). One observes 
that the eigenmodes of index r > 1 are all of first order in Knudsen number. In the asymptotic regime and to first 
order in the Knudsen number, J* 1 , a^ v on the r.h.s. of the relations (p~0|) can also be expressed in terms of any other 
eigenmode, e.g. the mth one, m > 1, such that for all n > 1, m ^ n, 

Y (i) 8 (i) 

^^^ + 0(Kn 2 ,KnR-\R- 2 ), 



Xn" fa 



y (2) 8 (2) 
n ~~ (2) R (2) m 



©(Kn^RnR^R^ 2 ) 



Then, choosing, e.g. m = 2 for £ = 1 and m = 1 for i — 2, in the above relations and using pjji 1 ' " w 
we can write, 



(11) 
(12) 

— Z^n=0 ' 



(2) ^iiv 



(2) 



fl (i) 

(1) Pn 
rn (!) 
,n=2 Xn 



A 2 — ii r0 A Q 



0(Kn) 



jv 2 



(2) Pn 



■rn ( 2 ) 
71=1 X" 



X{ v 



(2) 



O(Kn) 



(13) 



Using = n M and p$ v — , we replace the zeroth eigenmodes, X$ and Xq^, in Eq. (ITBl and express all irreducible 
moments, up to first order in Knudsen number, in terms of and , 



.(1) 



+ (kt - n^«o) ^2 = n#n" + 0(Kn) 
P2 



(2) 

"ro"^ + 2 ^2) - n%!m) xr = Wo>^ v + 0(Kn) , 



>(2) 



(14) 



where K r = X)fe=o #1 r rfc °4 ' ^ = J2k=o T rk a k > w ^ n r ^ = (-^ an d we se l ^00 = 1- This f° rmu l a is exact 
even for r = 0, for which the first-order terms in Knudsen number vanish. 

In Ref. [9(, X% and X£ u in Eqs. (fT4l are substituted by their Navier-Stokes values, i.e., 

B (1) B (2) 
X» = \l», xr = \a^. (15) 

x\ X\ 

These equations are then used to express Eqs. (0 in terms of only n M and (up to some order in Kn and R _1 ), 
leading to the following set of dynamical equations (in the massless limit) Q, 



M -J- = kP 1 + J» + /C7 + TV" , 

7T^ = 2^ + + KT + ft"" 



(16) 



In the above equations of motion, all nonlinear terms and couplings to other currents were collected in the tensors 
Jt*, K*, TV 1 , J» v , and W. The tensors and J"" contain all terms of first order in the product of Knudsen 
and inverse Reynolds numbers, 



(17) 



where we defined = V^, = V^P , and = (VV - V 1 '^) /2. The tensors /C' 1 and W v contain all terms 
of second order in Knudsen number, 



£" = Ki(7 A " / / 1 , + K2(T^ V F u + K 3 I"9 + RiF^O + KsU^Iv + Kg A^9jyCT + SyV^e, 

/C"" = r7 1 w A (A1 w l/ > A + r7 2 6lCT^+ 773^^^ +%<r< /i w , '> A + + fj 6 F^ F u > + fj 7 I^ F^ + ffc^ I v) + fj 9 V { "F u >. 

(18) 
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The tensors TZ^ and IZ^ contain all terms of second order in inverse Reynolds number, 

11" = ipin v TT^, TV" = (^TT^TT^ + iptfnP U "> . (19) 

These terms only appear when we consider nonlinear contributions in <5/k in the expansion of the single-particle 
distribution functions appearing in the collision term. Here, we will not consider such nonlinear contributions. In 
Eqs. (HH), terms of order C(Kn 3 ), 0(R _3 ), 0(Kn 2 R _1 ) and C(KnR~ 2 ) were omitted. 

Up to order C(R _1 ), O(Kn) and C(R_~ 1 Kn), Eqs. (fTo| have the same structure as those derived within IS theory 
[if all terms are kept in the derivation [HI]], although with different transport coefficients @: in principle, in RTRFD 
the transport coefficients carry information from all moments of the distribution function. In practice, the precision 
required for the values of the transport coefficients sets a limit on the number of moments that actually need to be 
taken into account. We found that, for the purposes of this work, 37 moments turn out to be sufficient. 

We now note that there is a crucial problem with Eqs. (1161) : the terms of higher order in Knudsen number, e.g. 
V</* J w > ) y<M p v ) ^ A£d v a av and V M 0, have second-order spatial derivatives and thus render these equations parabolic, 
i.e., acausal. Note that this happens despite the introduction of the slowest microscopic time scales as relaxation times 
for the particle diffusion current and the shear-stress tensor. Therefore, in order to solve Eqs. (|16[) in a relativistic 
setting, we have to neglect (at least some of) the terms (in) /C M , JC^ V . 

However, there are situations where these terms are important and should not be simply neglected. For instance, 
gradients of shear stress (the term ~ Rq A^9^cr Al/ C /C M ) and of the shear-stress tensor (the term ~ £ n „ Va7t a C J^) 
could be of the same order, resulting in source terms of similar magnitude in the equation for the particle diffusion 
current. Neglecting the former can spoil the agreement with numerical solutions of the Boltzmann equation, especially 
in the case where the Navier-Stokes term ~ 1^ is small compared to these second-order terms. Similar arguments 
apply to a situation where the gradients of 1^ and of the particle diffusion current are of the same order of magnitude 
and the Navier-Stokes term ~ a^ v is small. We shall study these situations in detail in Sec. IIVI 

We remark that this is also prone to happen in cases where the Navier-Stokes terms are of different order of 
magnitude. This would require the introduction of different Knudsen numbers for the various dissipative currents and 
the previous power-counting scheme has to be adapted to this situation. In this work, we do not pursue this avenue 
further from a formal point of view. Rather, we discuss in the next section how to include second-order terms in 
Knudsen number in a way that preserves hyperbolicity. 



III. CAUSAL TRANSIENT FLUID DYNAMICS UP TO SECOND ORDER IN KNUDSEN NUMBER 



The parabolic and, thus, acausal nature of the equations of motion (| L6[) can be understood as follows. The 
main assumption of RTRFD is to approximate the quickly varying eigenmodes of the Boltzmann equation by their 
asymptotic (i.e., Navier-Stokes) values. The approximation of the quickly varying eigenmodes with r > 1 happened 
in Eq. (|10l) . while the substitution of the eigenmodes with r = 1 (for 1=2) and r — 2 (for I — 1) by their 
Navier-Stokes values occurred in Eq. (|15p . However, in this step it was implicitly assumed that these eigenmodes 
relax instantaneously to their corresponding Navier-Stokes values, consequently leading to acausal behavior and, 
ultimately, to the parabolic terms in Eqs. (|16p . 

In order to obtain hyperbolic equations of motion which do not simply neglect terms of order 0(Kn 2 ), it is necessary 
to refrain from the substitution (|15l) . This can be simply done by keeping X% and in Eqs. (|14l) as independent 
dynamical variables instead of replacing them by their Navier-Stokes values. However, we do not solve differential 
equations of the type (HJ for these variables. Rather, for a given r we solve the first Eq. (fl4l) for Xl^ and, for a given s 
(not necessarily equal to r) the second Eq. (fT3j) for X£ v '. Thus, we can replace these variables by a set of irreducible 
moments p%, p% v . Without loss of generality, we choose r = 2 and s = 1, i.e., the irreducible moments p 1 ^ and p 1 ^" as 
independent dynamical variables. In this way, we obtain 

# = A^V + A$P2 + 0(Kn 2 , Kn RT 1 , R~ 2 ) , 

= AfoV* + A^Vr + 0(Kn 2 ,Kn R-\R~ 2 ) , (20) 

where we defined 

_ 0ffiw r -nffi«a ,(i) _ K n ( r ^ - K r , (2 ) _ »iqV - Offim ,(2) _ ^oVo ~Vr , 91 , 

"20 K — K 2 K0"20 _ K 2 "10 VO — Vl "10 ^0 ~ Vl 

These relations hold for any r and, therefore, up to order C(Kn 2 ,Kn R _1 ,R~ 2 ) it is possible to approximate all 
irreducible moments solely in terms of n/ 1 , P2, , and p±. 
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Note, however, that the relations (|20l) are only valid for the irreducible moments with positive r. Nev- 

ertheless, similar relations can also be obtained for the irreducible moments with negative r. This can be done by 
substituting Eq. ([3]) into Eq. (j4j , thereby expressing the irreducible moments with negative r in terms of those with 
positive r, 



rn rn ) 



(22) 



where we defined the following thermodynamic integral 



(2£- 



(23) 



Then, we approximate p^_ r (in terms of n M and p%) and p^ r (in terms of and p^ u ) by substituting Eqs. (|2H)) into 
Eq. (p2]l. leading to 



P-r 



LIU 

P-r = 



Ni 

V .F«A (1) 

/ y ^ rn ' v n0 



V .FWA (1) 

/ t J rn /x n2 
n=0,#l 



^ + 0(Kn 2 ,KnB7 1 ,B7 2 ), 



V- 77(2) x(2) 
/ , ■* rn ''nO 



N 2 



V- 77(2) A (2) 



+ 0(Kn 2 ,KnRr 1 ,R- 2 ). 



(24) 



Using these relations, Eqs. ([5]) can then be closed in terms of 7r MI/ , pj, and Pi"- In order to derive the equations 

of motion, we first multiply Eqs. (0 by Tnr and sum over r. Next, we use Eqs. ([20)) and (|24|) to replace all irreducible 
moments and p^ v appearing in the equations by n^, /?2 , Tr' 1 ", and p^ . The scalar irreducible moments and those 
with rank larger than two are replaced by zero. Additionally, all covariant time derivatives of ato, /3q, and are 
replaced by spatial gradients of fluid-dynamical variables using the conservation laws, d^N 11 = 0, d^T^ — 0. Then 
we obtain, 



f A M 



dn° 



dr 

dff a/3 



f ft u VIJj — 6 ri^O + l A* M/ d\w x —T if^F - A ft a^ v - A ir^I 



t^fP— + ir^ = 2ffa^ + 2t^oj^ x - 5^ v 6 - - f m f^ F v > + 4 n V<" ft v) + A ffn n^ J"> 



where we defined the vectors, 



P2 



Pi 



(25) 



(26) 



In this approximation, RTRFD becomes a theory with 21 moments as dynamical variables (had we included the scalar 
moments, there would have been 23 moments) while, in the previous approximation, i.e., Eqs. (fTB)) . there were only 
13 moments (including bulk viscous pressure, it would have been 14 moments). Above, f„, t^, £ n7r , £^ n , S nn , S wn , f n7r , 
fvrjr, f^n, A„ n , A nw , and X nn are 2x2 matrices, while k and rf are two-component vectors. The microscopic formulas 
for these transport coefficients, for a massless gas of particles, are shown in Appendix [3] Here, they are computed for 
a gas of classical particles with a constant cross section a, choosing AT X = 4 and N 2 = 3. The values for the diffusion 
and viscosity coefficients, k and if, and for the relaxation time matrices, f n and f n , are, 



Amfp'T-O 

1.295 -0.053/T 2 



0, 1596 \ y 
2.3616 T 2 ) ' A mfp P 



1.268 
6.929 T 



A 



mfp 



5.18 T z 



2.787 



A 



mfp 



0.912 0.136/T\ 
-3.647 T 2.456 ) 



(27) 



The transport coefficients of the nonlinear terms in the equation of motion for n M are 



&nn 
A m fp 



1.295 

5.18 T 2 



-0.0883/T 2 
4.645 



A,, 



0.524 -0.0341/T 2 \ A n7r _ / 0.1677/T 



Amfp 4P 



Amfp 

0.0973/T 2 
-2.6106 



2.096 T 



2.863 



p 



' A ; 



mfp 



A m f p 

-0.4723/T 0.0973/T 2 
13.111 T -2.611 



0.6708 T 



-0.0288/T 2 
-0.1147 



(28) 
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while those in the equation of motion for ir^ are 

K-k _4 / 0.912 0.17/T\ jw_ _ / 1.5688 0.2261/T\ jt™_ J_ / 0.2228 T 
A mfp " 3 ^ -3.647 T 3.0698 ) ' A mfp ~ ^ -6.2751 T 5.0956 ) ' A mfp ~ P V °- 8913 T * 

iwn_ _ ( 0.2228 T 0.0476/T\ _ / 0.1186 T 0.0084/T\ 

A mfp ~ ^ 0.8913 T 2 1.0096 ) ' A mfp ~ ^ -0.4744 T 2 -0.0338 J ' 

Note that these equations of motion are hyperbolic and include terms of C(Kn 2 ). We remark that, while this 
approach increases the domain of validity of the equations of motion without making them parabolic, it does not 
solve the intrinsic problem of the coarse-graining procedure: if one attempts to go to an even higher order in Knudsen 
number, the equations become once more parabolic. This can be solved in a similar fashion, by again increasing the 
number of irreducible moments describing the system. 



0.0714/T \ 
1.5144 ) ' 

(29) 



IV. COMPARISON WITH THE BOLTZMANN EQUATION 



In order to test the validity of RTRFD, we compare it to the solution of the Boltzmann equation, similarly as in Ref. 
Q. The numerical method of choice for solving the Boltzmann equation is the Boltzmann approach for multiparton 
scattering (BAMPS) [ll|, while the macroscopic field equations are solved using the viscous SHarp And Smooth 
Transport Algorithm (vSHASTA) [l4j . Both fluid dynamics and the Boltzmann equation are solved in Cartesian 
coordinates with a flat space-time metric = diag(l, —1, —1, —1). 

We consider a (3+l)-dimensional system, but assume that the matter is homogeneous in the y — z plane, al- 
lowing for an inhomogeneous matter distribution only in the longitudinal x-direction. This effectively leads to a 
(1+1) -dimensional problem. We consider two different initial conditions. In case I, the system is initialized with a 
homogeneous fugacity, Ao = e a ° = 1, but with an inhomogeneous pressure profile in the longitudinal direction. In 
practice we smoothly connect two temperature states T(_ 



0.4 GeV and T, 



(+c 



0.25 GeV via the Woods-Saxon 

parametrization with thickness parameter D = 0.3 fm. In case II, the pressure Pq = gT^^" 1 /n 2 is homogeneous and 
the fugacity distribution is also given by a Woods-Saxon profile with D = 0.3 fm, interpolating between A(_oo) = 1 and 
A(+oo) — 0.2. For the degeneracy factor we use g = 16. In both cases, matter is initialized in local thermodynamical 
equilibrium, i.e., with all dissipative currents (and eigenmodes of the Boltzmann equation) set to zero, and at rest, 
i.e., with a vanishing collective velocity — 0. These initial conditions arc shown in Fig.Q] 
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FIG. 1. (Color online) Initial conditions for cases I and II. 



In both cases, we consider two exemplary values for the cross section a = 2 and 8 mb, and consider the solutions 
after the system has evolved for 6 fm in time. We compare the solution of the Boltzmann equation with that of 
traditional IS theory [including terms omitted in the original work Q but quoted in Ref. [HI], as well as with RTRFD 
at various levels of approximation. Equations (|16[) contain 13 moments as independent dynamical variables (and 14, 
if we include bulk viscous pressure). The calculation of the transport coefficients in these equations can be done with 
increasing accuracy, as more irreducible moments are considered in the expansion ([3]). The lowest possible accuracy 
is reached if no more than the original 13 (14) irreducible moments are considered for the calculation of the transport 
coefficients. At the next level, we include one more set of irreducible moments of tensor-rank one and two (and one 
more scalar moment in the case of non- vanishing bulk viscous pressure), which leads to a total of 21 (23) irreducible 
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moments. In this way, the number of irreducible moments entering the transport coefficients increases by 8 (9) at each 
successive level of approximation. For the purpose of this work, we found that going to the third level of iteration, 
i.e., considering 13 + 8 • 3 = 37 moments (14 + 9 • 3 = 41 in the case of nonvanishing bulk viscous pressure) is sufficient 
to reach the desired accuracy in the values of the transport coefficients. In the following, we shall compare RTRFD 
with 13 dynamical degrees of freedom and with the transport coefficients computed with 13 and with 37 moments. 
We shall term these variants of RTRFD "13/13" and "13/37". In addition, we also solve Eqs. ([25]). These contain 
21 dynamical degrees of freedom. We compute the corresponding transport coefficients using 37 moments. We shall 
refer to this variant of RTRFD as "21/37". 

In the following figures, the numerical solutions of the Boltzmann equation shall always be displayed by open dots, 
the results of IS theory by black dash-dotted lines, the solution of RTRFD "13/13" by green dashed lines, that of 
RTRFD "13/37" by blue dotted lines, and that of RTRFD "21/37" by solid red curves. 

In order to verify the different fluid-dynamical theories discussed in this paper the solutions of the Boltzmann 
equation must be calculated to a very high precisison. For this purpose we performed 5 • 10 4 BAMPS runs and 
computed the fluid-dynamical quantities as averages over these runs. 

In Fig. [2] we show the fugacity (top) and thermodynamic pressure (bottom) and in Fig. [3] the heat flow = 

— (e + Po)n^ /n (top) and shear-stress tensor (bottom) for case I. The Boltzmann equation and the fluid-dynamical 
theories were solved for a = 2 mb (shown in the left panels of each figure) and for a — 8 mb (shown in the right 
panels). For a = 8 mb, the thermodynamic pressure and shear-stress tensor computed in all fluid-dynamical theories 
are in good agreement with the numerical solutions of the Boltzmann equation. As we decrease the cross section we 
expect the agreement between macroscopic and microscopic theory to become worse. This explains why, for a = 2 
mb, the pressure and shear-stress tensor computed within fluid-dynamical theories deviate more strongly from those 
computed via the microscopic theory. Nevertheless, compared to the fugacity and heat-flow profiles, the agreement is 
not too bad, even for the smaller value of the cross section. 

The initial pressure gradient in case I drives, via conservation of momentum, the creation of velocity gradients. On 
the other hand, the gradient of fugacity is initially zero and turns out to remain small throughout the evolution. In 
this situation, higher-order terms involving gradients of velocity and of the shear-stress tensor, e.g. ReA^d^a^ C /C M 
and ^n^-A^V ' C in the particle diffusion equation (TTB| . become of the same order as the respective (first- 
order) Navier-Stokes term kI^. Therefore, if terms of this type are not properly taken into account, we expect large 
deviations from the solution of the Boltzmann equation. This can be seen in Figs. [5] and [3] when comparing IS theory, 
RTRFD "13/13", as well as RTRFD "13/37" with the Boltzmann result. In all of these variants, the parabolic 
term ~ k 6 is either absent (IS theory and RTRFD "13/13") or has to be dismissed (RTRFD "13/37") for reasons 
of causality. In addition, IS theory and RTRFD "13/13" do not have the correct value for £ n „, because we did not 
include a sufficiently large number of irreducible moments in its computation. Although RTRFD "13/37" features 
(within the desired accuracy) the correct value for this transport coefficient (as well as for R e ), it does even worse in 
describing the fugacity and heat-flow profiles than the previous two theories. This is because the term ~ R e could not 
be taken into account for reasons of causality, although it is of the same order of magnitude as the term ~ £ n7r . These 
problems of fluid-dynamical theories with only 13 dynamical variables are resolved by RTRFD "21/37" which is the 
only fluid-dynamical theory considered here that contains all contributions of second-order in the Knudsen number 
in a hyperbolic fashion. 

In Fig. [4] we show the fugacity (top) and thermodynamic pressure (bottom) and in Fig. [5] the heat-flow (top) and 
shear-stress tensor profiles (bottom) for case II. As before, the Boltzmann equation and the fluid-dynamical theories 
considered were solved for a — 2 mb (shown in the left panels) and for a — 8 mb (shown in the right panels). Again, 
we expect, and see, better agreement between fluid dynamics and the Boltzmann equation for the larger value of the 
cross section. While the fugacity profiles are in good agreement with the solution of the Boltzmann equation for all 
fluid-dynamical theories and both values of the cross section, the heat flow is not well described in IS theory and in 
RTRFD "13/13": IS theory predicts values for the heat flow which are smaller in magnitude than the Boltzmann 
equation, while RTRFD "13/13" predicts larger values, even for a = 8 mb. On the other hand, both RTRFD "13/37" 
and RTRFD "21/37" describe the heat flow very well or even perfectly, respectively, for both values of the cross 
section. The reason is that the diffusion coefficient k has the correct value in these theories (while it deviates by 

- 30% in both IS theory and RTRFD "13/13"). 

Since, in case II, the initial pressure gradient is zero and turns out to remain small throughout the evolution, the 
velocity gradients remain small as well. In this situation, it is important to include higher-order terms that couple the 
shear-stress tensor to heat flow. This is the reason why the solutions of IS theory and RTRFD "13/13" (where these 
higher-order terms vanish in the massless limit) are not in good agreement with that of the Boltzmann equation for the 
thermodynamic pressure and the shear-stress tensor, for both values of the cross section. On the other hand, RTRFD 
"13/37" does a better job in matching the Boltzmann equation. It is not perfect, because the higher-order terms ~ fjs 
and ~ fjs were dropped. The best agreement is, again, found within RTRFD "21/37" where all second-order terms 
in the Knudsen number are taken into account. 
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FIG. 2. (Color online) Fugacity and thermodynamic pressure profiles at t = 6 fm for case I, for a = 2 mb (left panels) and 
cr = 8 mb (right panels) . 



Note that, in Fig. [3] the BAMPS results for the shear-stress tensor are strongly fluctuating. This happens because, 
in this special case, the values of the shear-stress tensor are of the same magnitude as the statistical fluctuations in 
BAMPS. In order to reduce the statistical fluctuations and to achieve a better resolution, a significantly larger amount 
of runs would be required. 



V. SUMMARY AND CONCLUSION 



In this paper, we compared the equations of motion of RTRFD at various levels of approximation with numerical 
solutions of the Boltzmann equation for two different types of initial conditons (labeled case I and II). Also, we showed 
how to use the RTRFD formalism to derive equations of motion that are hyperbolic and, at the same time, include 
terms up to second order in the Knudsen number. By a careful comparison with numerical solutions of the microscopic 
theory, we demonstrated that this formalism is able to handle problems with strong initial gradients in pressure or 
particle number density. 

The initial conditions, cases I and II, were chosen in such a way that considerably different spatial gradient profiles 
are generated throughout the fluid-dynamical evolution. In case I, the pressure gradient is initially large which gives 
rise to large velocity gradients in the later stages of the evolution. This means that, in case I, the shear-stress 
tensor is mainly generated by its corresponding Navier-Stokes term, i.e., by gradients of velocity. On the other hand, 
the fugacity gradient is initially zero in case I, and remains relatively small throughout the evolution of the fluid. 
Therefore, in case I the heat flow is not mainly created by its Navier-Stokes term, i.e., by the gradient of fugacity, 
but by the coupling term to the shear tensor and shear-stress tensor, i.e., the terms A^V\tt^ and A^V^c^ m Eqs. 
(|16p . Therefore, in this case the higher-order terms in Knudsen number must be included and one really needs to 
solve the hyperbolic equations derived in this paper to obtain a good agreement. The fact that IS theory is always 
deviating from the microscopic theory when it concerns heat flow means that it does not predict correctly the terms 
of order one and two in Knudsen number. On the other hand, in RTRFD "13/37" and RTRFD "21/37" all transport 
coefficients are computed with a sufficiently large number of irreducible moments. This guarantees that all terms of 
the desired order are included and is the reason for the better agreement of these fluid-dynamical theories with the 
microscopic theory. The reason why RTRFD "13/37" fails in certain situations is that important terms have to be 
neglected in order to preserve hyperbolicity and causality. 

In case II, the fugacity gradient is initially large while the pressure gradient is zero. This means that the heat flow 
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originates mainly from its Navier- Stokes term while the shear-tress tensor originates mainly from its coupling to heat 
flow, i.e., the terms V <Al n" > , V <A '/ l/> , I <f *I" > , and n <f *I u> in Eqs. HT5J). The fact that the heat flow calculated 
from IS theory deviates from the solution given by the microscopic theory even in this case, is evidence that the 
Navier-Stokes term of this theory does not contain the correct transport coefficient. The coupling of the shear-stress 
tensor with the heat flow in IS theory is also not correctly taken into account. 

In conclusion, the resummation of irreducible moments for the computation of the transport coefficients was essential 
to obtain a good agreement with the microscopic theory. It provides not only the correct values for the shear viscosity 
and heat conduction coefficients, but also for the transport coefficients that couple the respective dissipative currents. 
Moreover, in situations where higher-order terms are important, one has to make sure to include them in a hyperbolic 
way, and not simply drop relevant contributions because they are parabolic. These two factors resolved the previously 
observed differences between the solution of IS theory and of the Boltzmann equation observed in Ref. Q . 

As expected, and explicitly demonstrated in this paper, the agreement between solutions of RTRFD and the 
Boltzmann equation depends on the value of a. For the cases considered in this paper, we obtained a good agreement 
for (7 — 8 mb, while for a = 2 mb we started to notice small deviations. In order to improve the agreement for smaller 
values of the cross section, we would have to include more moments of the Boltzmann equation to describe the state 
of the system, i.e., such moments would have to contribute not only to the values of the transport coefficients but 
also as independent dynamical variables. We leave this investigation to future work. 
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FIG. 4. (Color online) Shear-stress tensor and heat flow at t = 6 fm for case II, for a = 2 mb (left panels) and a = 8 mb (right 
panels) . 

Appendix A: Transport coefficients 

In this appendix we list all transport coefficients of fluid dynamics calculated in this paper. The microscopic 
formulas for the diffusion and viscosity coefficients, k and if, and for the relaxation time matrices, f n and f T , are 





The transport coefficients of the nonlinear terms in the equation of motion for ft 1 * are 

Ni 

3 



1 Nl 



r=0,^l 



3r (1) A (1) 



3t (1) A (1) 



5r (1) A (1) 



5r (1) A (1) 




(-> r w (1) A (2) 
{l - r) T Qr A r _ 10 

{l - r) T 2r A r _ 10 



(1)\(2) 



(2-r)r^A 



(2-r)r^A 



(1)\(2) 



(1)^(2). (2) 

T 00 ■'lr A r l 
(1)W2) X (2) 

'20 lr A rl 



(1)^(2) ,(2) 
'00 •'lr A r0 

(1)^(2). (2) 
T 20 -'lr A r0 



(1)^(2). (2) 
T 00 •'lr A rl 

(l) 7 -(2) x (2) 
T 20 •'lr A rl 



4P 



Ni 



Jl)r 
T 0r ir+2,1 

r=0#l V T 2r 




r (D A (2) 
'Or A r-l,l 
r (D A (2) 



(Al) 



(A2) 
(A3) 
(A4) 



(A5) 



12 




while those in the equation of motion for n 1 * 1 ' are 
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